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We study numerically the outcome of the phase separation instability of a dual-species Bose- 
Einstein condensate in an optical lattice. When only one excitation mode is unstable a bound pair of 
a bright and dark soliton-like structures periodically appears and disappears, whereas for more than 
one unstable mode a persistent soliton-antisoliton pair develops. The oscillating soliton represents a 
regime where the two-species condensate neither remains phase-separated nor is dynamically stable. 
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A two-component atomic Bose-Einstein condensate 
(BEC) is a veritable laboratory of dynamics. Experi- 
mentally observed long lifetimes of two-species conden- 
sate mixtures have prompted many studies of phase sep- 
aration 0- 0, 0, HI, , symmetry-breaking transitions 
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and shock £ [waves, and vector soliton struc- 
12l , 13, 3|. Similar phenomena are also 



H, spin 
tures [Gj, _ 

possible in nonlinear optics I15|. A two-species BEC may 
be loaded in an optical lattice [l6| , which should afford an 
even wider tunability and control of system parameters 
than non-lattice systems. 

Phase separation is a generic phenomenon that oc- 
curs in different forms of matter and plays an impor- 
tant role, e.g., in mixing of liquids. It has previously 



been shown [9|, [ill, |l3j, [14|, [15| that phase separation be 



tween two superfluid components in a BEC can stabilize 
a bright soliton that would otherwise not occur. In the 
phase separation regime a bright solitary wave in one 
component exists trapped inside a density dip of a co- 
propagating dark soliton in the other component. Gener- 
alizations of dark-bright solitons, such as 3D particle-like 
solitons in a two-species BEC lfl, 12, 17 1 and analogous 
cosmic strings in the early Universe cosmology (Til ], are 
similarly stabilized by phase separation. 

In this Letter we study theoretically the outcomes of 
a dynamical (modulational) phase separation instability 
of a two-component BEC in a lattice using the number 
of unstable excitation modes as the organizing principle. 
For one unstable mode, we surprisingly find a quasiperi- 
odically emerging bound pair of bright and dark soliton- 
like structures (a soliton and an antisoliton) as the first 
sign of the phase separation instability (Fig. [1]), whereas 
more than one unstable mode can lead to a formation of a 
persistent vector soliton or to more familiar violent phase 
mixing dynamics. The quasipcriodic soliton formation 
represents a regime where the system neither is stable, 
nor can undergo phase separation. It is reminiscent of 
a pulsating instability due to anomalous dispersion in a 
single-species BEC in a lattice 



19, 2fJ . However, here the 




FIG. 1: (Color online) Creation of a pair of a bright soliton 
(peak) and a dark soliton (dip) resulting from a dynamical 
instability. The parameters are xn = 1, X22 = 100, X12 = 12, 
and Ji = J2 =0.5, with the lattice size L — 64. The density 
of the dark soliton has been rescaled to make it visible. The 
arrows indicate the directions of lattice sites (n), time (t) and 
the density of each species (\i/)n^ | 2 )- 



not occur in a single BEC. We are led to speculate that 
an instability unfolding as pulsations is a widespread phe- 
nomenon in nonlinear systems with one unstable mode. 
Our proposed lattice system could be used as a testing 
laboratory for the pulsating instability and soliton forma- 
tion; it is of particular interest that current experiments 
should be able to demonstrate these phenomena. 

In order to study in detail the dynamics of a two- 
species BEC in a regime where the uniform density 
profile in the lattice becomes dynamically unstable we 
performed linear stability analysis of a classical mean- 
field dynamics of the BECs along the lines of Ref . 21 1 . 
A weakly interacting BEC in a deep optical lattice 



can be described by the tight-binding version of the 
Gross-Pitaevskii equation, known as discrete nonlinear 
Schrodinger equation (DNLSE). The two-species version 
of DNLSE reads 



dt Vn 



CO 



fc=l 



nature of the instability is phase separation, which can- 



where Jj (Jj 



> 0) denotes the hopping amplitude of 
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the atoms between adjacent lattice sites for the species 
j = 1, 2, and the interaction coefficients Xjk describe the 
onsite atom-atom interactions that are proportional to 
the intra-species and the inter-species s-wave scattering 
lengths djj and ajk (j ^ k), as well as to the over- 
lap of the lowest vibrational state wave functions 
(the Wannier functions) of the two species in each lat- 
tice site, so that Xjj — iirhNiCLjj /m i Jd 3 r\4>H ) \ 4 and 
Xj k - 2nh^Nj^a jk /^Jd 3 r\^ ) \ 2 \ ( j ) i n k) \ 2 , for j £ k. 
Here nij is the atomic mass of the j th BEC component 
and fi = m\m2j (jn\ + tti-i) is the reduced mass. We 
use the normalization Yln=i IV'n ^ | 2 = 1 as the number of 
atoms is conserved for each species separately. Through- 
out our calculation we consider a finite lattice length L, 
with periodic boundary conditions such that the lattice 
site with the index n = is equivalent to the lattice site 
with the index L. Technically we are dealing with a ring 
lattice, but our results clearly should have counterparts 
for other types of boundary conditions. 

The DNLSE flU may be solved with the standard 
ansatz i/}„ — e 1 ^™-^ 4 ) /y/L t providing a station- 
ary solution with the two chemical potentials u>j = 
— 2 Jj cos pj + Y^t=iXjk/L. Here pj is the carrier wave 
momentum for the j th species. The boundary conditions 
quantize the pj as pj = Pj(2w/L), where Pj is an inte- 
ger that may be chosen from the interval [— L/2, L/2). 
For notational convenience we take the number of lattice 
sites L to be even. 

The stability of the steady state solution for the two- 
species BEC may be studied using the classical Bogoli- 
ubov expansion. By substituting 

= (l + u tf>e<9" - [ u (J')]* e - t ' ? ™)e l( ^ n -^ t) /v / I (2) 

into Eq. (TT]) and collecting the first-order terms in Uq 

(i) 

and Vq , which are assumedly small, we obtain a matrix 
equation of the type [2l[, id£/dt = 9Jt£, where 971 is a 4 x 
4 non-Hermitian matrix and £ = [uq 1 \vq 1 \u < ^\v q 2 * > ] T . 
The eigenvalues of 9JI correspond to the normal mode 
frequencies of the system. Due to the periodic boundary 
conditions the quasimomcnta q are also quantized as 



flow, ti 



'1,9 



q and 



q = 2nQ/L , 



(3) 



with Q taking on the integer values except zero in the 
range [—-§,-§)■ A simple form for the analytic expres- 
sions for excitation frequencies Q q are obtained when 
two BECs have the same atomic currents, Jisin(pi) = 
J 2 sin(p 2 ) [II, namely, 



0„ 



ti ± 



(4) 



where t\ = 2 J\ sin(pi) sin(g) represents a Dopplcr shift 
term of the excitation frequencies due to the supcrfluid 



*3 = \/( V X,q - V l,q) 2 + 16ei, <3 e2,gX? 2 / i2cOS Pl C0S P2 , 

(5) 

arc defined in terms the single-condensate normal mode 
frequencies Vj q (without the Doppler shift term) 



3,9 ~ 0,1 



C0S (Pj)[^j,g cos(p 3 ) + 2xjj/L] , (6) 



and tj. q = 4Jj sin 2 (g/2) is the spectrum of an ideal, 
non-moving BEC. The flow is stable if the frequencies 
in Eq. (j4|) are real for all 5 / 0, otherwise there are 
small excitations in the system that grow exponentially 
in time. For the purposes of the stability analysis, the 
excitation modes with the indices Q and —Q may, and 
will, be regarded as equivalent. 

For a uniform single-component BEC with repulsive 
atom-atom interactions the condensate flow in a lattice 
turns dynamically unstable when the flow velocity ex- 
ceeds a critical value 22j. In a two-component BEC dy- 
namical instabilities may arise both from critical velocity 
and from phase separation of the two species [2l|. Here 
we study the dynamics in the normal dispersion regime, 
for quasimomcnta pi,P2 < § and interaction parameters 
Xjk > 0, so that the instability is ascribed exclusively 
to phase separation. According to Eq. ([4]), those lowest 
energy modes turn unstable for which q satisfies 

X12 > XiiX22+q 2 L(j 2 Xii +J1X22) /2+q i L 2 j 1 j 2 /4:. (7) 

Here we are studying the lowest excitation modes with 
< Q <C L, where q is given by Eq. (O, and write 
jk = Jk cos pk- The number of unstable modes may be 
obtained from Eqs. and by solving for Q. In 
the limit of a large number of lattice sites L — > 00, we 



obtain the condition for instability X12 > X11X22 [21|. 



This is an obvious counterpart of the iconic condition for 
the instability of free-space condensates with equal-mass 
atoms 0, H[, <xi2 > ana22- A similar analysis may also 
be performed in the case that one of the BECs exhibits 
anomalous dispersion due to superfluid flow, but the re- 



sulting stability conditions are significantly altered [21 



We next simulate the dual-species dynamics in a lattice 
using a sixth-order accurate split-step FFT method as 
described in [23[ for a single-condensate Gross-Pitaevskii 
equation. We start with a homogeneous non-moving mix- 
ture of two BECs, adding a small amount of random noise 
(~ 10~ 6 ) to the condensate amplitudes to speed up the 
instability. We first concentrate on what we call weak 
instability, when only one excitation mode with Q = 1 is 
unstable, and also assume that the interaction strengths 
satisfy X22 3> Xn- Because of the proliferation of param- 
eters, we do not wish to use any of them in particular as 
a scale of dimensional quantities. Instead, we have arbi- 
trary compatible units for energy, frequency, and time, 
so that dimensional quantities appear as dimensionlcss 
numbers in our further analysis. 
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In Fig. [T] wc show the time evolution of the den- 
sity profiles \ipn | 2 of initially non-moving BECs as the 
instability runs its course. The parameter values are 
Xi = 1) Xi = 100 and X12 = 12, for L = 64 and 
J\ = J-2 = 0.5. The critical value of inter-species in- 
teractions for the onset of the instability for these pa- 
rameters is X12 = 10.75, according to Eq. (JT]). The value 
X12 = 12 in Fig. [1] corresponds to one unstable mode, 
as the threshold for the second mode to become unsta- 
ble is X12 = 12.75. Here one expects phase separation, 
and a region with a pronounced perturbation between 
the densities (numbers of atoms per lattice site) of the 
two species indeed develops. However, Fig. Q] shows that 
upon further evolution the species mix again. The system 
returns close to its initial state, subsequently an abortive 
phase separation occurs again, and so on. 

In this feature the system bears a strong resemblance 
to the previously demonstrated pulsating instability that 
is associated with superfluid flow [3, ■ The mecha- 
nism is also likely to be as we surmised before: At the 
point of instability the stationary state, which used to 
be an elliptic fixed point of the dynamics, turns into a 
hyperbolic fixed point. In the process a homoclinic or- 
bit emerges that starts from the hyperbolic fixed point 
in the unstable direction and returns to the same fixed 
point in the direction of a corresponding excitation mode 
that shrinks exponentially. Starting in the vicinity of the 
hyperbolic fixed point, the system approximately tracks 
the homoclinic orbit, which leads to periodic recurrences 
of phase separation and the initial state. 

It may seem counterintuitive that the minority species, 
the one with much weaker interspecies interactions, may 
form a local density maximum in spite of the repulsive 
interactions and normal dispersion, but there is a well- 
known analog in the literature in the form of a bound pair 
of bright and dark solitons both in free space [HI and in 
a harmonic trap 
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131 ]. In what is ultimately a case of 
phase separation, the background of the coincident dark 
soliton stabilizes the bright soliton, providing an effective 
trap due to the repulsion between the two species when 
the variation in the total atom density is energetically less 
favorable than in relative atom densities [H, E3- The 
pulsating variations of atom density we have described 
thus far neither constitute an energy minimum nor are 
they stationary, but on the strength of the analogy we 
nevertheless refer to our finding as a bound pair of a 
dark and a bright soliton. 

The large ratio X22/X1I! e -g-j d ue to population im- 
balance, is essential for the emergence of solitons. For 
nearly equal nonlinearities we do not observe recogniz- 
able localized structure, and the structure narrows with 
the increasing ratio X22/Xn- O n the other hand, hav- 
ing more than one excitation mode dynamically unsta- 
ble eventually leads to persistent phase separation, for 
X22/X11 S> 1 in the form of a dark-bright soliton, with 
the associated phase kink in the dark soliton. We illus- 
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FIG. 2: (Color online) Top row from left: Profiles of persis- 
tent soliton pairs for (xn, X22, X12) = (1,10,8), (1,100,20), 
and (1,1000,60) for the common parameters J\ — J2 — 0.5, 
L = 64. Here each example starts with four unstable modes. 
Bottom row from left: Persistent soliton pairs for X12 = 15 
(two unstable modes), 36 (eight unstable modes), and 20 
(eight unstable modes). Here we have J\ = J2 = 0.5, xn — 1, 
X22 = 100, and L — 64, except in the last case L = 256. The 
dip in the antisoliton (dashed line) has been multiplied by a 
factor of 2.5 to enhance the visibility. The x axis in each panel 
counts the lattice sites whereas the y axis is for the densities 
of minority (full line) and majority (dashed line) species. 



trate this in Fig. [5] by showing the density profiles of 
the final soliton pair for three different ratios of X22/X11 
(top row). In each case the number of initially unsta- 
ble modes was four. The persistent soliton profiles also 
become narrower when the number of unstable modes 
increases; Fig. [2] (bottom row). A large number of unsta- 
ble modes can result in the two species being mixed in a 
more complicated way with no recognizable solitons 2J] , 
analogously to the previous experiments 0]. Moreover, 
we found that, in the case of several unstable modes, in- 
creasing the number of sites leads to more solitons, as 
also demonstrated in Fig. O For balanced populations 
we may also observe maximally phase-separated density 
profiles where the two species occupy the separate halves 
of the lattice [25| . All soliton profiles in Fig. [5] are sur- 
rounded by sound wave excitations appearing as a fluc- 
tuating background density. While these soliton pairs 
appear stable, there is nothing in our translationally in- 
variant model to pin them, so they jitter around. 

This dark-bright soliton composite can be viewed as 
a bound state of two particles, a vector soliton, as is 
seen by setting it in motion. Here we resort to numeri- 
cal experiments inspired by oscillatory soliton dynamics 
that has been both discussed theoretically [ll| and ob- 
served experimentally Initially we have a homoge- 
neous plane wave with pj = for each component. Wc 
switch on a weak trapping potential in the axial direc- 
tion, V n = 7}n(n — -§-) 2 with k = 10~ 6 , and wait for the 
pulsating soliton pair to emerge. For the parameter val- 
ues we have used it always comes at the trap center, and 
stays there if the trap remains on unchanged. However, 
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FIG. 3: Motion of bright (top) and dark (bottom) pulsating 
solitons in an added harmonic potential V n = \n(n — -§-) 2, 
Lattice position from a total of 64 sites (vertical) is shown 
as a function of time (horizontal). Except for the harmonic 
trap with k = 1CP 6 , the parameters are the same as in Fig.Q] 
Brighter shading represents higher densities in the top panel, 
and lower densities at the bottom. 

once we have the soliton pair, we abruptly displace the 
center of the trap, whereupon the trapping forces make 
the solitons oscillate back and forth. For the pulsating 
soliton pair, the oscillation amplitude in the trap can be 
much larger than the initial displacement. This is demon- 
strated in Fig. [3] where the trap center was displaced by 
one lattice site, yet we see oscillation amplitudes of about 
20 sites. The bright and dark solitons in a pulsating pair 
stay together while they oscillate back and forth just like 
the components of a persistent vector soliton, but for a 
persistent soliton the oscillation amplitudes are compa- 
rable to the initial displacement. 

In the case of a single species, weak instability requires 
a weak nonlincarity. The pulsating instability in a single- 
component supcrfluid flow therefore entails long time 
scales, and the attendant experimental problems. In 
contrast, the present pulsating soliton- antisoliton pair 
seems to be within easy reach of current experimental 
techniques. This is because the two species stabilize 
each other, and the interactions can be large while the 
system remains in the regime of weak instability with 
just one unstable excitation mode. The time scale of 
the instability may thereby be shortened and, e.g., is 
r = l/|Syfi g | ~ 80 for the parameters of Fig. [TJ If the 
root energy scale is the recoil energy, as is reasonable for 
optical lattices, for 87 Rb this translates to 20 ms. Con- 
tinuing with our qualitative example, the time between 
two consecutive pulses in Fig. [T] would be 0.4 s. 

One could also think of in-principle problems such as 
failure of the DNLSE to approximate well the underlying 
quantum field theory. However, the DNLSE only fixes 
the nonlincaritics \ki oc ^Nk\fW~i and, for a given values 
of xm the classical mean-field theory in an optical lattice 
becomes increasingly accurate for higher atom numbers. 
In our example estimates of quantum effects on the pul- 
sating instability in a single-species BEC we have found 
that the DNLSE may fare well for atom numbers as small 
as hundreds [2fj| . 

In conclusion, we have described theoretically a fea- 
ture that uniquely combines phase separation, soliton- 



antisoliton pair, and pulsating instability: In the limit of 
only one unstable excitation mode a bound pair of bright 
and dark solitons that emerges and disappears nearly pe- 
riodically while a dynamical instability of a two-species 
BEC in an optical lattice runs its course. Basically, we 
have repeated attempts at phase separation that do not 
quite make it. However, if more than one excitation mode 
is unstable, phase separation goes through and produces 
a persistent soliton-antisoliton pair. The dual-species 
condensate in an optical lattice also exhibits other insta- 
bilities besides phase separation, and additional peculiar 
behaviors arc likely when these instabilities take off. 
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